;;=====================================================================================
;; calculate RH pdf 
;;=====================================================================================
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl" 
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl" 

begin

 plotname = "out"
 plotform = "pdf"

 PDF1D = True
 PDF2D = False ;; True 
 PDF1D = False ;;True 
 PDF2D = True

 if(PDF1D) then
    plotname = plotname + "_1D"
 end if

 if(PDF2D) then
    plotname = plotname + "_2D"
 end if

 StringFontHeightF  = 0.035
 cnLineLabelFontHeightF = 0.025
 tmXBLabelFontHeightF = 0.030 ;;25
 tmYLLabelFontHeightF = 0.030 ;;25
 tiMainFontHeightF  = 0.030
 tiXAxisFontHeightF = 0.035 ;;0
 tiYAxisFontHeightF = 0.035 ;;0
 lbLabelFontHeightF = 0.025

 fna = "mod_sparticus_wsub_ni_pdf_25.nc" 
 fnb = "mod_sparticus_wsub_ni_pdf_25.nc" 

 fla = addfile(fna,"r") 
 flb = addfile(fnb,"r") 

 pdf1_0 = fla->pdf1_0 ;; Ni < 0.01  #/L outside cloud  

 pdf1_0 = where(pdf1_0.eq.0.,-999,pdf1_0)

 pdf1_0@_FillValue = -999 ;;where(pdf2_1.eq.0.,-999,0) 

 pdf2_0 = flb->pdf2_0 ;; Ni < 0.01  #/L outside cloud  
 
 pdf2_0 = where(pdf2_0.eq.0.,-999,pdf2_0) 

 pdf2_0@_FillValue = -999 ;;where(pdf2_1.eq.0.,-999,0) 

;*************************************************************
; plots
;*************************************************************

  plot = new (6,"graphic")

  wks  = gsn_open_wks(plotform,plotname)


  if(PDF1D) then 

  res  = True
  res@gsnDraw                = False
  res@gsnFrame               = False
  res@xyLineThicknessF       = 2

  ;; string 
  res@gsnLeftStringFontHeightF   = StringFontHeightF
  res@gsnCenterStringFontHeightF = StringFontHeightF
  res@gsnRightStringFontHeightF  = StringFontHeightF

  ;; tickmark
  res@tmXBLabelFontHeightF = tmXBLabelFontHeightF
  res@tmYLLabelFontHeightF = tmYLLabelFontHeightF

  res@tmXBMode = "Explicit"	
  res@tmXBValues = (/0,3,6,9,12,15,18,21,25/) 
  res@tmXBLabels = (/"0","3","6","9","12","15","18","21","25"/) 

  ;; axis
  res@tiXAxisFontHeightF = tiXAxisFontHeightF
  res@tiYAxisFontHeightF = tiYAxisFontHeightF

  res@gsnLeftString        = "a) Ni < 0.01 # L~S~-1~N~"
 ;res@gsnCenterString      = ">0.01 #/L"
  res@tiYAxisString        = "Frequency (%)"
  res@tiXAxisString        = "" ;;"RHice (%)"
  plot(0) = gsn_csm_xy (wks, pdf1_0@bin_center, pdf1_0, res) 

  resP                     = True                ; modify the panel plot
  resP@txString            = ""
  ;;resP@gsnPanelRowSpec = True                   ; tell panel what order to plt
  gsn_panel(wks,plot,(/3,2/),resP)


  end if 


  if(PDF2D) then 

  colormap = "WhiteBlue"
  colormap = "wh-bl-gr-ye-re"
  colormap = "BlueYellowRed"

  gsn_define_colormap(wks,colormap)

  res    = True
  res@gsnDraw              = False    ; do not draw picture
  res@gsnFrame             = False    ; do not advance frame
  res@cnFillOn             = True     ; turn on color fill
  res@cnFillMode           = "RasterFill"       ; Raster Mode
  res@cnLinesOn            = False    ; no contour lines
  res@cnLineLabelsOn       = False    ; no contour line labels
  res@cnInfoLabelOn        = False
  res@gsnSpreadColors      = True     ; use full colormap
 ;res@cnLabelBarEndStyle   = "ExcludeOuterBoxes"
  res@lbOrientation        = "vertical"
  res@lbLabelAutoStride    = True
  res@tiXAxisString        = "Temperature (K)" 
  res@tiYAxisString        = "wsub (cm s~S~-1~N~)"
 ;res@tiXAxisOffsetYF      = -0.05 ;;0.105

  ;; string 
  res@gsnLeftStringFontHeightF   = StringFontHeightF
  res@gsnCenterStringFontHeightF = StringFontHeightF
  res@gsnRightStringFontHeightF  = StringFontHeightF

  ;; tickmark
  res@tmXBLabelFontHeightF = tmXBLabelFontHeightF
  res@tmYLLabelFontHeightF = tmYLLabelFontHeightF

  ;; axis
  res@tiXAxisFontHeightF = tiXAxisFontHeightF
  res@tiYAxisFontHeightF = tiYAxisFontHeightF


  res@cnLevelSelectionMode =  "ManualLevels"
  res@cnMinLevelValF       = 0.
  res@cnMaxLevelValF       = 0.2 
  res@cnLevelSpacingF      = 0.02 
  res@cnMaxLevelValF       = 0.6
  res@cnLevelSpacingF      = 0.02
  res@cnMaxLevelValF       = 0.5 
  res@cnLevelSpacingF      = 0.05;;1

  ;;res@tmXBMode = "Explicit"
  ;;res@tmXBValues = (/ 1.e-3,2.e-3,5.e-3,1.e-2,2.e-2,5.e-2,1.e-1,2.e-1,5.e-1/) 
  ;;res@tmXBLabels = (/ 1.e-3,2.e-3,5.e-3,1.e-2,2.e-2,5.e-2,1.e-1,2.e-1,5.e-1/) 

  res@tmXBMode = "Explicit"
  res@tmXBValues = (/ 1.e-1,1.,1.e1,1.e2,1.e3,1.e4,1.e5,1.e6/)
  res@tmXBLabels = (/ 1.e-1,1.,1.e1,1.e2,1.e3,1.e4,1.e5,1.e6/)


  ;res@lbLabelFontHeightF  = 0.010               ; make labels smaller


  res@gsnLeftString        = ""
 ;res@gsnCenterString      = ">0.01 #/L"
  res@tiXAxisString        = "Ni (# L~S~-1~N~)" 
  res@tiYAxisString        = "wsub (cm s~S~-1~N~)"
  plot(0) = gsn_csm_contour (wks, pdf2_0, res)

  resP                     = True                ; modify the panel plot
  resP@txString            = "" 
 ;;resP@gsnPanelRowSpec = True                   ; tell panel what order to plt
  gsn_panel(wks,plot,(/2,2/),resP)



  end if 


;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; 
;; 
;; scatter plot 
;; 
;;
;; plot = new(4,graphic)  
;;
;; wks  = gsn_open_wks("pdf","RHaim_inside")              ; specifies a ps plot
;; 
;; res                     = True                   ; plot mods desired
;; res@gsnDraw              = False        ; do not draw picture
;; res@gsnFrame             = False        ; do not advance frame
;; res@gsnMaximize         = True                   ; maximize plot in frame
;; res@xyMarkLineModes     = (/"Markers"/)  ; choose which have markers
;; res@xyMarkers           = 16                     ; choose type of marker 
;; res@xyMarkerColor       = "red"                  ; Marker color
;; res@xyMarkerSizeF       = 0.005                  ; Marker size (default 0.01)
;; res@xyDashPatterns      = 1                      ; solid line 
;; res@xyLineThicknesses   = (/1/)                ; set second line to 2
;;
;; res@tiMainString        = "SPARTICUS - outside cloud"  ; title
;; res@tiMainString        = "SPARTICUS - inside cloud"  ; title
;;
;; plot(0)  = gsn_csm_xy (wks,px,py,res)        ; create plot
;; plot(1)  = gsn_csm_xy (wks,px,pa,res)        ; create plot
;; plot(2)  = gsn_csm_xy (wks,px,pb,res)        ; create plot
;; plot(3)  = gsn_csm_xy (wks,px,pc,res)        ; create plot
;;
;; ;; draw panel with white space added
;; resP                 = True
;; resP@gsnPanelYWhiteSpacePercent = 5
;; resP@gsnPanelXWhiteSpacePercent = 5
;; gsn_panel(wks,plot,(/2,2/),resP)
;;
;;  
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;; 

 
end






